/*
 *  lib-smv: Matrix - Vector Library targeted to Structural Mechanics Problems
 *  Copyright (C) 2006-2008  St.Brcic, Lj.Zugic
 *
 *  This program is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  This program is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU General Public License for more details.
 *
 *  You should have received a copy of the GNU General Public License
 *  along with this program.  If not, see <http://www.gnu.org/licenses/>.
 */
 
/** @file EigVal.h
 *  @brief Definition of class EigVal - solution of the eigenvalue problem
 *
 *  EIGENVALUE PROBLEMS:  Standard and Generalized   
 *  A*x = L*x and A*x = L*B*x
 *
 *  Matrices A and B are real, square symmetric and positive definite (SPD).
 *
 *  Implemented routines are:
 *  (1) Jacobi transformation (only for standard problems)
 *  (2) Householder transformation combined with QL transformation
 *  (3) Lanczos algorithm (for standard and generalized problems) 
 */ 

#ifndef _EIGVAL_H
#define _EIGVAL_H

#include <iostream>
#include <fstream>
#include <iomanip>
#include <cstdlib>
#include <cstdio>
#include <cmath>
#include <cassert>
#include "Arrays.h"
#include "iArrays.h"
#include "Jacobi.h"
#include "TQLI.h"
#include "Lanczos.h"
#include "Cholesky.h" 

using namespace std;

namespace smv_lib
{
//! usefull enumerations for eigenvalue analysis
//! Complete or Partial solution
enum EV_ANAL   {COMPLETE, PARTIAL};
//! Standard or Generalized problem
enum EV_TYPE   {STANDARD, GENERAL};
//! Standard Jacobi, Householder or Lanczos method
enum EV_STANDARD {sJACOBI, HOUSE, sLANCZOS};
//! Equivalent Standard, Generalized Jacobi or Lanczos method
enum EV_GENERAL {eqSTAND, gJACOBI, gLANCZOS};

class EigVal
{
private:
	//--- basic input data
	//! Matrices A and B (square real and positive definite)
	Matrix A, B;
	//! Order of matrices A / B
	int N;
	//! Number of eigenvalues and/or eigenvectors
	int EV_Number;

	//--- solution parameters
	//! Complete or Partial solution
	EV_ANAL EV_Anal;
	//! Standard or Generalized problem
	EV_TYPE EV_Type;
	//! Standars Jacobi, Householder or Lanczos method
	EV_STANDARD EV_Standard;
	//! Equivalent Standard, Generalized Jacobi or Lanczos method
	EV_GENERAL EV_General;

	//--- control data
	//! Indicator: true if matrix A is given
	bool isMatrix_A;
	//! Indicator: true if matrix B is given
	bool isMatrix_B;
	//! Indicator: true if number of eigenvalue/vector pairs is given
	bool isNumber;
    //! Indicator: true if Standard problem is considered
    bool isStandard;
    //! Indicator: true if Generalized problem is considered
    bool isGeneral;
    //! Indicator: true if method of solution is given
    bool isMethod;
	//! Indicator: true if the name of the eigenvalue problem is given
	bool isNameGiven;
	//! Indicator: true if the soltion is obtained
	bool isSolution;

	//--- basic output data
	//! Modal matrix (columns of X are eigenvectors)
	Matrix X;
	//! Array of Eigenvectors (i.e. each xi is a column 'i' of matrix X)
	Vector* xi;
	//! Vector of Eigenvalues (lambda)
	Vector Lam;
	// Vector kruznih frekvencija (stampa se u datoteku radi kontrole)--30.10.07.
    Vector Omega;
   	// Vector perioda oscilovanja (stampa se u datoteku radi kontrole)--30.10.07.
    Vector T;
    // Vector frekvencijaa (stampa se u datoteku radi kontrole)--30.10.07.
    Vector f;

	//! Name of the eigenvalue problem (output file is Name + ".txt")
	string Name;

	//! solution methods for the eigenvalue problem
	//! Jacobi transformation
	Jacobi*  m_Jacobi;
	//! Householder transformation and QL transformation
	TQLI*    m_TQLI;
	//! Lanczos algorithm
	Lanczos* m_Lanczos;
    //! Symmetrix matrix - Cholesky decomposition---30.10.07.
    Cholesky* m_Chol;
    //! Indicator: true if Cholesky LLT is selected---30.10.07.
    bool isCholesky;
	

public:
	//! empty (default) constructor 
	EigVal(void);
	//! constructor with ONE matrix and a number of eigenvalue/vector pairs
	EigVal(Matrix& a, int n = M_DEFAULT);
	//! constructor with TWO matrices and a number of eigenvalue/vector pairs
	EigVal(Matrix& a, Matrix& b, int n = M_DEFAULT);
	~EigVal(void);
	//-------------------- set and get functions (in case of a default constructor)
	//! set Matrix A (for standard problem)
	void setMatrix_A(Matrix& a);
	//! set Matrices A and B (for generalized problem)
	void setMatrix_B(Matrix& a, Matrix& b);
	//! set number of eigenvalue/vector pairs
	void setNumber(int num);
	//! set control parameters
	void setAnal(EV_ANAL anal)      { EV_Anal = anal; }
	void setType(EV_TYPE type)      { EV_Type = type; }
	void setStand(EV_STANDARD st)   { EV_Standard = st; isStandard = true; }
	void setGeneral(EV_GENERAL gen) { EV_General = gen; isGeneral = true; }
	//! set the name of the eigenvalue problem
	void setName(string name);

	//! get number of eigenvalue/vector pairs
	int getNumber() const          { return EV_Number; }
	//! get control parameters
	EV_ANAL getGroup() const       { return EV_Anal; }
	EV_TYPE getType() const        { return EV_Type; }
	EV_STANDARD getStand() const   { return EV_Standard; }
	EV_GENERAL getGeneral() const  { return EV_General; }

    //! solution of the problem (main routine)
	void solveEV_Problem(string fileName = "");
	//---- solution results
	//! get the modal matrix (eigenvectors)
	Matrix getEV_Matrix();
	//! get all eigenvalues
	Vector getEV_Spectar();
	//! get eigenvector number i
	Vector getEV_Vector(int i);
	//! get the eigenvalue number i
	double getEV_Value(int i);
    //! write results to an output file
    void writeEV_Results(string fileName = "");
private:
	bool checkInput();
	void selectMethod();
    // implemented algorithms for standard problem
    void standardJacobi();
    void Householder_TQLI();
    void standardLanczos();
    // implemented algorithms for generalized problem
    void equiStandard();
    void generalJacobi();
    void generalLanczos();
	//void elmhes (double **a, long n);
	//void hqr    (double **a, long n, double wr[], double wi[]);
};

}   // end of namespace smv_lib

#endif //  _EIGVAL_H

//=============================== end of file
